#Required Packages
library(sf)
## Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1
library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v ggplot2 3.3.5 v purrr 0.3.4
## v tibble 3.1.4 v dplyr 1.0.7
## v tidyr 1.1.3 v stringr 1.4.0
## v readr 2.0.1 v forcats 0.5.1
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(sp)
library(spdep)
## Loading required package: spData
library(raster)
##
## Attaching package: 'raster'
## The following object is masked from 'package:dplyr':
##
## select
## The following object is masked from 'package:tidyr':
##
## extract
library(leaflet)
library(RColorBrewer)
#Prepare Data
#Lab 2
#County Shapefile
counties <- sf::read_sf("../data/CBW/County_Boundaries.shp") %>% sf::st_make_valid()
#BMP data
bmps <- read_csv("../data/CBW/BMPreport2016_landbmps.csv")
## Rows: 69601 Columns: 18
## -- Column specification --------------------------------------------------------
## Delimiter: ","
## chr (11): StateAbbreviation, GeographyName, Geography, Agency, BMPShortName,...
## dbl (7): AmountSubmitted, AmountBackedOut, AmountNotBackedOut, AmountCredit...
##
## i Use `spec()` to retrieve the full column specification for this data.
## i Specify the column types or set `show_col_types = FALSE` to quiet this message.
#Lab 3
#filtered States
states.filtered <- sf::read_sf("../data/states_filter.shp")
#Lab 4
#recent little wabash river shp file
lwr_new <- sf::read_sf("../data/LWR_Line.shp") %>% sf::st_make_valid()
#1940 little wabash river shp file
lwr_1940 <- sf::read_sf("../data/Left 1940.shp") %>% sf::st_make_valid()
##counties in IL shp file
IL_co <- sf::read_sf("../data/IL_BNDY_County_Py.shp") %>% sf::st_make_valid()
#DEM little wabash river
lwr_ras <- raster::raster("../data/lwr_ras.tif")
#Task 1
#total cost of BMPs funded by county
bmps_trim <- bmps %>% mutate(., FIPS.trimmed = stringr::str_sub(GeographyName, 3, 5)) %>% group_by(FIPS.trimmed) %>% summarise(tcost = sum(Cost, na.rm = T))
#Join with county shp
counties_bmps <- left_join(counties, bmps_trim, by = c("COUNTYFP10" = "FIPS.trimmed"))
counties_bmps$tcost <- as.integer(counties_bmps$tcost)
#set basemap variable
Base_Map <- providers$CartoDB.Positron
#find equal interval breaks
c.max <- max(counties_bmps$tcost, na.rm = T)
c.min <- min(counties_bmps$tcost, na.rm = T)
br <- (c.max - c.min)/5
start <- 0
bins <- c()
for (i in 1:4) {
start <- start+br
bins <- append(bins, start)
}
#set bins and pal variables
pal <- colorBin(c("#F9E79F", "#7DCEA0", "#17A589", "#2E86C1", "#154360"), domain = counties_bmps$tcost, bins = c(0, bins, Inf), na.color = NA)
#map frame
base <- leaflet(counties_bmps) %>%
setView(lng = -76.55906218576634, lat = 40.507740442861356, zoom = 6) %>% addProviderTiles(Base_Map, group = "Base_Map")
#Hover Label
labels <- sprintf(
"<strong>%s</strong><br/>$ %s total cost of BMPs",
counties_bmps$NAME10, counties_bmps$tcost
) %>% lapply(htmltools::HTML)
#choropleth map
base %>%
addPolygons(
fillColor = ~pal(tcost),
weight = 2,
opacity = 1,
color = "white",
fillOpacity = 0.85,
label = labels,
labelOptions = labelOptions(
direction = "bottom",
style = list(
"color" = "green",
"font-family" = "serif",
"font-style" = "italic",
"box-shadow" = "3px 3px rgba(0,0,0,0.25)",
"font-size" = "12px",
"border-color" = "rgba(0,0,0,0.5)"
)
)
) %>%
addLegend(
pal = pal,
values = ~tcost,
title = "Total Cost of BMPs by County",
position = "bottomright"
)